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Abstract 

Recently, the paradigm that the dynamic magnetosphere displays sandpile- 
type phenomenology has been advanced, in which energy dissipation is by 
means of avalanches which do not have an intrinsic scale. This may in turn 
imply that the system is in a self organised critical (SOC) state. Indicators 
of internal processes are consistent with this, examples are the power law 
dependence of the power spectrum of auroral indices, and in-situ mag- 
netic field observations in the earth's geotail. An apparent paradox is 
that, rather than power laws, substorm statistics exhibit probability dis- 
tributions with characteristic scales. Here we discuss a simple sandpile 
model which yields for energy discharges due to internal reorganization 
a probability distribution that is a power law, whereas systemwide dis- 
charges (flow of "sand" out of the system) form a distinct group whose 
probability distribution has a well defined mean. We analyse the model 
over a wide dynamic range whereupon two regimes having different in- 
verse power law statistics emerge, corresponding to reconfigurations over 
two distinct scaling regions: short scale sizes sensitive to the discrete na- 
ture of the sandpile model, and long scale sized up to the system size 
which correspond to the continuous limit of the model. The latter are 
anticipated to correspond to large scale systems such as the magneto- 
sphere. Since the energy inflow may be highly variable, we examined the 
response of the model under strong or variable loading and established 
that the power law signature of the large scale internal events persists. 
The interval distribution of these events is also discussed. 



1 Introduction 

Recently there has been interest in relating observed characteristics of global 
energy transport in space plasmas to "sandpile" models which dissipate energy 
by means of avalanches (Consolini, 1997; Chapman et ai, 1999). When such 
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models exhibit scale free, inverse power law statistics in the probability distri- 
butions of energy released by avalanches, and of avalanche length and duration, 
they are candidates for description in terms of self organised criticality (SOC) 
(Bak et al., 1987, 1988; Lu, 1995) see also Jensen (1998) and references therein. 
The power spectra may also have an inverse power law ("l/f) signature, and 
SOC was introduced to explain the ubiquity of such spectra and of fractality in 
Nature. 

Chang's suggestion (Chang, 1992, 1998a, 1998b) that the magnetosphere 

is in an SOC state has motivated application of avalanche models to the so- 
lar wind-magnetosphere-ionosphere system. Observational motivation includes 
sporadic nature of energy release events within the magnetotail ("Bursty Bulk 
Flow" events (Angclopoulos et al, 1996)), power law in-situ magnetic field power 
spectra (Hoshino et al., 1994), and power law features of magnetospheric index 
data, notably AE which is an indicator of energy dissipated by the magneto- 
sphere into the ionosphere. Tsurutani et. al. (1990) described a broken power 
law AE spectrum; indicative of SOC but not conclusive as power law power 
spectra are not unique to SOC systems (Jensen, 1998) . Consolini (1999) has 
recently used AE data taken over a ten year period to construct the distribu- 
tion D(s) of a burst measure s extending the result obtained for one year in 
(Consolini, 1997). This work strongly suggests that inverse power law burst 
statistics are a robust feature of the AE data, albeit with an exponential tail 
and some evidence of an additional lognormal component. It is currently less 
clear that the power law is solely of intrinsic magnetospheric origin and it may in 
fact be related to the behaviour of the solar wind energy input (Freeman et al., 
1999). Since the global reconfigurations of the magnetotail (substorm events) 
appear to have occurrence statistics with a well defined mean. Chapman et al. 
(1998) demonstrated a simple avalanche model (another example is (Pinho and 
Andrade, 1998)) that in principle has relevance for the magnetosphere in that 
it yields systemwide avalanches where the statistics have a well defined mean 
(intrinsic scale), whereas their internal avalanche statistics are scale free. 

An additional consideration for space plasmas are the means by which the 
conjecture of scale free statistics can be tested. We wish to test the hypothesis 
that the probability distributions of energy dissipated, length scales and dura- 
tion of avalanches are of power law form as seen in the original observations 
of SOC (Jensen, 1998) in a slowly driven sandpile. Since to test for power law 
statistics we need to maximise the range of sizes of events, the required statis- 
tical experimental evidence requires long runs of data. In the magnetospheric 
system this will then imply that both the instantaneous value and smoothed 
local mean of the loading rate (the solar wind) will have strong variation. 

In this paper we illustrate some robust features of the avalanche statistics of 
the simple avalanche model (Chapman et al., 1998, 1999) which are needed for 
application of such a model to space plasma data. We investigate to what extent 
the model gives inverse power law avalanche statistics under slow loading. We 
show how these statistics are modified under strong and/or variable loading. We 
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shall also see that two distinct regimes for energy transport, both with power 
law avalanche statistics but with different slopes, emerge from the sandpile 
algorithm, depending on the size of the system. In addition we show how these 
may be characterized in terms of the interval distribution of events on different 
spatial scales. 

2 A two- regime avalanche model 
2.1 The algorithm 

Sandpile algorithms generally include an array of nodes, at each of which there 
is a variable amount (height) of sand; a critical gradient (difference in height 
between neighbouring nodes) which, if exceeded by the actual gradient, triggers 
local redistribution of sand; and algorithms for redistribution and fuelling. The 
main measured output is the statistics of the emergent avalanche distribution. 
Kadanoff et al. (1989) gave an early classification of such models. The relation- 
ship of the models to experimental sandpiles, and to the ideal concept of SOC 
remain topics of active research (see for example Jensen (1998), also Dendy and 
Helander (1997)). 

The sandpile model used here is described in more detail in (Chapman et 
al., 1999; Hclandor et al, (1999)). We have a one- dimensional grid of N equally 
spaced cells one unit apart, each with sand at height hj and local gradient 
Zj = hj — hj+i. There is a repose gradient zr below which the sandpile is always 
stable, and with respect to which heights hj and the gradients Zj arc measured. 
Each cell is assigned a critical gradient Zcj- If the local gradient Zj exceeds 
this, the sand is redistributed to neighbouring cells and iteration produces an 
avalanche. The critical gradients on each of the N nodes are selected from 
a uniform ("top hat") probability distribution P{zcj). P{zcj) is generated by 
choosing the Zcj at random with uniform probability from the range [a, h], and 
the integral of P{zcj) over all Zcj is unity. 

Sand is added at cell 1 at a rate g, and we normalize length and time to 
the mean loading rate. If the critical gradient is exceeded at cell 1, the sand is 
redistributed. Our redistribution rule is conservative and instantaneous: sand 
will propagate to cell 2 and if the local critical gradient there is exceeded, to 
cell 3 and so on. Within this avalanche, the sand is instantaneously "flattened" 
back to the angle of repose at which the sandpile is always stable. Propagation 
of an ongoing avalanche from one cell (fc) to the next (/c + 1) thus occurs if 



hk - hk+i > Zck 



(1) 



In consequence a quantity of sand A is deposited on the next cell: 



K+i = hk+i + A 



(2) 
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so that the gradient at k relaxes to the angle of repose (here normalized to zero) 

K - K+i = ZR = o (3) 

Since all cells within the ongoing avalanche 1,2..., k arc at the angle of repose 
following this conservative redistribution of sand to the k + 1 cell, we require 
the heights of all these cells to become: 

K..k = hi..k - J (4) 

The iterative procedure described above (where superscript * denotes an inter- 
mediate step) continues until the avalanche reaches a cell where the gradient is 
subcritical. The critical gradients at cells within the flattened post-avalanche 
region are then rerandomized as above. More sand is added then at cell 1 
until this end cell again becomes unstable, triggering another avalanche. An 
avalanche may be entirely an internal rearrangement of sand or may continue 
until it spreads across all A'' cells. We call the latter a systemwide discharge, in 
which case the entire sandpile is emptied and returns to the angle of repose. 

This relaxation rule allows propagation of information (correlation) across 
the avalanche and therefore potentially on all length scales in the sandpile. This 
permits the possibility of scale free self organising behaviour in a one dimensional 
system, whereas other SOC reorganisation rules (Bak et al., 1987; Jensen, 1998) 
have typically required two dimensions. 

The total energy dissipated by an avalanche (internal or systemwide) is just 
the difference in the potential energy in the entire sandpile before and after the 
avalanche 

JV N 
C^e = ^ ^1 \before " X/ \after (5) 

i=i i=i 

An example of a time series for the energy is shown in Figure 1. The 5000 
cell sandpile was loaded slowly {g = 0.001) with respect to the mean value of the 
Zcj ■ The critical gradients were uniformly and randomly distributed in the range 
[0.5, 1.5]. With the angle of repose normalized to zero, the time evolution is by 
systematic growth as sand is added, interspersed with systemwide avalanches 
where the energy falls back to zero, and internal avalanches where the energy is 
reduced to some nonzero value. 

The statistics of the energy released in internal and systemwide avalanches 
for two longer runs of this sandpile under more realistic conditions of fluctuating 
input are shown in Figure 2, with the normalized probability distribution P{de) 
for both internal and systemwide avalanches plotted as a single population. As 
in all sandpile runs here, the populations comprise over 5 x 10^ internal and 
2 X 10"^ systemwide avalanches. 

The two runs in figure 2 are for slow (< g >= 0.001, diamonds) and fast 
(< g >= 10, circles) mean loading rates, giving an indication of the expected 
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behaviour of a system with strong variation in the driver, such as the solar wind 
driven magnetosphere (see also (Watkins et al., (1999)). For both values of 
inflow rate internal avalanches show distinct inverse power law regions with a 
turndown at small de, whereas the systemwide events which have a characteris- 
tic mean (Chapman et al, 1998, 1999) appear as a bump at the high energy end. 
The distinct behaviour of the systemwide avalanches, independent of the inflow 
rate, shown here is a necessary condition for applicability to the magnetosphere 
(Chapman et al, (1998)). The internal avalanches show different behaviour un- 
der slow and fast loading. Essentially, the smaller events are destroyed as we 
increase the average loading rate, making larger events more probable. Hence 
the normalized probability of larger events shows an increase on the plot. Im- 
portantly, their power law slope (here of index —1) is preserved and is thus a 
robust feature that should be apparent in observations under variable loading. 
We now discuss the internal avalanches in more detail. 

2.2 System scales and power law index 

Two lines a{de)~^ and 0{de)~^ are drawn on Figure 2 (and all subsequent 
figures). The values a = 0.25, 7 = 0.65 and = 2 are an approximate fit through 
the points. Under slow loading the sandpile exhibits two distinct regimes. In the 
case where P{zcj) = S{z — a), (a any positive constant) the sandpile evolution 
with time has been found analytically and, if the system is normalized to have 
total length unity, it can be shown that P(de) = {de)~^ (Helander et al., 1999). 
The arguments of Helander et at, (1999)) lead us to expect a region of power 
law index —1 in P(fie) for sandpile with P{zcj) of finite width. We might also 
anticipate that as the width of P{Zcj) is decreased more of the total range of 
P{de) would be characterized by a power law index —1, but surprisingly this is 
not so. 

The P{de) for four sandpile runs are superposed in Figure 3, differing only by 
the choice of P{zcj)- For three of these the same mean critical slope < Zcj >= 1 
but three different widths (0.01,0.1,1) have been used. The fourth also has 
P{zcj) with width 0.1, but has a different mean < Zcj >= 0.1. In this latter 
case we have rescaled de ^ dex 100 since, on average, the heights of sand needed 
for instability will be smaller by an order of magnitude, so that (5) will yield 
values of de that are on average smaller by two orders of magnitude. Figure 3 
suggests that all features of the probability distribution are robust against the 
choice of P{zcj), which effectively represents the local condition for instability. 

Avalanches dissipating smaller amounts of energy might be expected to ex- 
tend over smaller length scales, illustrated in Figures 4-6 where we replot the 
data shown in Figure 3, showing only the contribution from successively longer 
avalanches. Avalanches with lengths > 1,8,64 respectively are shown. 

Independent of the details of P{zcj) we see that the power law index ~ —0.65 
corresponds to avalanches that extend over less than ~ 64 cells. Figure 4 then 
shows that the drop at de ~ 1 in Figure 3 corresponds to avalanches that are 
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one cell in length. 

The sandpile thus has three distinct regimes in its statistics: single cell 

avalanches that (as expected) are not power law; avalanches smaller than ^ 64 
cells, with power law index ~ —0.65, which may reflect the discrete nature of 
the grid; and avalanches longer than ~ 64 cells and up to the system size, with 
power law index —1, which may approach a continuous limit for the system. 

The interval distribution (that is, the probability distribution of time in- 
tervals between events) provides further insight into the difference between 
avalanches on the small and the large scale. In figures 7-10 wc show a scries of 
interval distributions from the single sandpile run with constant slow fuelling 
g = 0.001 and with top hat probability distribution for the critical gradients 
Zc = [0.95, 1.05]. The figures show the time intervals between sucessively larger 
avalanches, that is, of all avalanches of length > 1,4,8 and 64. Since sand is 
always added at cell 1, and hence instability always occurs first at cell 1 a plot 
of the interval distribution for all avalanches (not shown) simply reproduces the 
probability distribution for Zc, ie, the the distribution of times between sucessive 
avalanches, all of which are triggered at the first cell. As we selectively plot the 
distributions of time intervals between longer avalanches, that have propagated 
further down the sandpile, we see the effect of the interaction of more cells in 
the sandpile. 

Figure 7 shows intervals between all avalanches that have propagated beyond 
cell 1 (ie length > 1). Here we see two characteristic timescales corresponding 
to avalanches that stop at cell 2, and those propagating beyond cell 2. Time 
is normalised to the inflow rate (such that unit sand is added to the sandpile 
in unit time) so that avalanches that reach cell 2 and stop will only occur after 
sufficient sand has been added to exceed the critical gradient at cell 1, which 
has mean value 1. Hence the minimiim time interval in this case is > 1. As we 
increase the miminum avalanche length considered, the minimum time interval 
also increases correspondingly. The detailed behaviour becomes complex for 
lengths jl but less than 64, that is, including avalanches which dissipate energy 
according to power law index ^ —0.65. The general trend however is for an 
increasing number of characteristic time intervals to appear as wc only consider 
avalanches of increasing length. When we exclude avalanches that dissipate en- 
ergy according to power law index ~ —0.65 by only considering avalanches of 
length > 64 (figure 10) the interval distribution has become continuous with 
cutoff at time interval 65 as we would expect. The large scale avalanches identi- 
fied as those dissipating energy with probability distribution that is power law 
index —1 therefore correspond to this continuous limit. Unlike laboratory plas- 
mas (see Chapman et ai, (1999)) these large scale events are expected to be 
relevant to astrophysical systems and are expected to be the robust observable 
in the case of the magnetosphere which has strong, variable loading. 
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3 Conclusions 



A simple one dimensional sandpile model has been developed with two distinct 
characteristics in the probability distribution of energy discharges. For inter- 
nal reorganisation there are two distinct inverse power law regimes, whilst for 
systemwide discharges (flow of sand out of the system) the probability distribu- 
tion has a sharply- defined mean. Our model may be applied to magnetospheric 
dynamics (Chapman et al., (1998)), for example in reconciling the apparent 
paradox of power law indexes in internal dynamics with substorm event statis- 
tics which have peaked distributions. 

Under slow loading the internal dynamics exhibits two regimes which have 
inverse power law statistics of index ~ —0.65 and —1, corresponding to re- 
configurations on distinct length scales. Short length scales may arise from the 
discrete nature of the grid, while we also see longer scales, up to the system size, 
that effectively approach a continuous limit of the model. We find a transition 
between these regimes at avalanche lengths of about 64 cells. 

For space plasma systems observations taken over long periods are required 
to test for possible inverse power law statistics. The loading of the system (in 
the case of the magnetosphere, the solar wind) is often characterised by both 
strong variability about a mean, and a large dynamic range of mean energy 
input. The inverse power law form of the statistics of large internal avalanches 
has been shown to be robust under fast loading. The effect of large loading 
rates is to exclude events which dissipate; small amounts of energy, which in 
our model results in a single inverse power law regime with downturn at lower 
energies. We expect such inverse power law avalanche distributions to be a 
persistent feature in long runs of data that include "fast" inflow conditions if 
the underlying system is governed by SOC. 
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Figure 1: The time evolution of energy in the 5000 ceh sandpile, for fuehing 
g = 0.001 and probabiUty distribution for the critical gradients that is top hat 
in the range [0.5, 1.5]. 
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Figure 2: The probability density of all internal and systemwide avalanches for 
a 5000 cell sandpile with variable fueUing < g >= 0.001 (o) and < g >= 10 (o) 
, with probability dis tribution for the critical gradients that is top hat in the 
range [0.5, 1.5]. 
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Figure 3: The probability density of all internal avalanches for a 5000 cell 
sandpile with constant fuelling g = 0.001 and four different runs with prob- 
ability distributions for the critical gradients that are top hat: o = [0.5,1.5], 
* = [0.95,1.05], o = [0.995,1.005] and (with rescaling de de x 100) 
+ = [0.05,0.15]. 
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Figure 4: The probability density of internal avalanches of length greater than 1 
for a 5000 cell sandpile with constant fuelling g = 0.001 and four different runs 
with probability distributions for the critical gradients that are top hat: o = 
[0.5, 1.5], * = [0.95, 1.05], o = [0.995, 1.005] and (with rescaling de ^ de x 100) 
+ = [0.05,0.15]. 
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Figure 5: The probability density of internal avalanches of length greater than 
8 for a 5000 cell sandpile with constant fuelling g = 0.001 and four different 
runs with probability distributions for the critical gradients as in the previous 
figure. 
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Figure 6: The probability density of internal avalanches of length greater than 
64 for a 5000 cell sandpile with constant fuelling g = 0.001 and four different 
runs with probability distributions as in the previous figure. 



15 



500 cells P(time between internal avalanches l>1) vz t 

0.1 2 I , 1 , , , , , 




1.5 2 2.5 3 3.5 4 4.5 5 5.5 6 6.5 

dt 



Figure 7: The probability distribution of time intervals between avalanches of 
length ^1, with time normalized to the constant fuelling rate. The fuelling 
g = 0.001 and probability distribution for the critical gradients that is top hat 
in the range [0.5, 1.5]. 
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Figure 8: The probability distribution of time intervals between avalanches of 
length ^4, with time normalized to the constant fuelling rate. The fuelling 
g = 0.001 and probability distribution for the critical gradients that is top hat 
in the range [0.5, 1.5]. 
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Figure 9: The probability distribution of time intervals between avalanches of 
length ^8, with time normalized to the constant fuelling rate. The fuelling 
g = 0.001 and probability distribution for the critical gradients that is top hat 
in the range [0.5, 1.5]. 
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Figure 10: The probability distribution of time intervals between avalanches 
of length ^64, with time normalized to the constant fuelling rate. The fuelling 
g = 0.001 and probability distribution for the critical gradients that is top hat 
in the range [0.5, 1.5]. 
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